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Thermo—inertial bouncing of a relativistic collapsing sphere: A numerical model 
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We present a numerical model of a collapsing radiating sphere, whose boundary surface undergoes 
bouncing due to a decreasing of its inertial mass density (and, as expected from the equivalence 
principle, also of the “gravitational” force term) produced by the “inertial” term of the transport 
equation. This model exhibits for the hrst time the consequences of such an effect, and shows 
that under physically reasonable conditions this decreasing of the gravitational term in the dynamic 
equation may be large enough as to revert the collapse and produce a bouncing of the boundary 
surface of the sphere. 

PACS numbers: 04.40.-b; 04.40.Dg; 95.30.Lz; 97.60.-s. 


I. INTRODUCTION 

In the study of gravitational collapse of massive stars, the inclusion of dissipative processes (in particular neutrino 
emission) is enforced by the fact that they provide the only plausible mechanism to carry away the bulk of binding 
energy, leading to a neutron star or black hole Q]. On the other hand, in cores of densities about 10^^g cm~^ the 
mean free path of neutrinos becomes small enough as to justify the use of diffusion approximation Hi This seems to 
be confirmed by the observational data collected from supernova 1987A, which indicates that the radiation transport 
regime prevailing during the emission process, is closer to the diffusion approximation than to the streaming out limit 

i _ 

Motivated by the comments above, in a recent paper Q , the Misner and Sharp approach to the study of adiabatic 
gravitational collapse @ was extended as to include dissipation in, both, the streaming out and diffusion approximation 
(for the case of pure free streaming approximation see Q)- Then from the comling of the dynamical equation to 
a causal transport equation in the context of Miiller-Israel-Stewart theory [a @1 it was obtained that the effective 
inertial mass density of a fluid element and the gravitational force term in the dynamical equation, reduce by a factor 
which depends on dissipative variables. This reduction, in its turn, might lead to the bouncing of the collapsing 
sphere, as discussed in p. 

As can be seen from inspection of the transport equation, such an effect is directly related to the presence of the 
inertial term Ta/s in the transport equation. This explains why we refer to such a bouncing as “thermo-inertial”. 

It is our purpose in this work to present a numerical model of a radiating collapsing sphere, where the above 
mentioned effect produces the bouncing of the boundary surface of the sphere, for physically acceptable values of all 
variables. 

Since we are mainly concerned with time scales of the order of magnitude of (or even smaller than) the hydrostatic 
time scale, as in the quick collapse phase preceding neutron star formation, we cannot rely on the quasistatic approx¬ 
imation, and therefore the full dynamic description has to be used [iaii3. This implies that we have to appeal to a 
hyperbolic theory of dissipation. The use of a hyperbolic theory of dissi pat ion is further justified by the necessity of 
overcoming the difficulties inherent to parabolic theories (see references |l2l | -[ ^ and references therein). 

The plan of the paper is as follows. In the next section we define the conventions and present the dynamical 
equation coupled to the transport equation. The model to be considered as well as the strategy for the numerical 
integration is presented in Section III. Finally, a discussion of results is presented in Section IV. 
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II. THE DYNAMICAL EQUATION OF THE DISSIPATIVE FLUID 

We consider a spherically symmetric distribution of collapsing fluid (for simplicity we shall consider the pressure to 
be locally isotropic) undergoing dissipation in the form of heat flow, bounded by a spherical surface S. We assume 
the interior metric to S to be comoving, shear free for simplicity, and spherically symmetric, accordingly it may be 
written as 


ds 2 ^ -A^{t, r)dt'^ + sin^ 

and hence we have for the four velocity U“ and the heat flux vector q°‘ 

U“ = 9 “ = . 

Then it can be shown that the following equation can be found from Bianchi identities 


+ P)DtU = -(^ + P) [m + AttPR^] - E'^DrP - E 


+ BDtq 

IX 


( 1 ) 

( 2 ) 

( 3 ) 


where is the energy density, P the pressure. 
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j__a 

Am' 


( 4 ) 


the proper radial derivative Dr, constructed from the radius of a spherical surface, as measured from its perimeter 
inside S, being 


with 


D =J-A 

^ R'dr' 


R = rB, 


( 5 ) 

( 6 ) 


and where dots and primes denote derivatives with respect to t and r respectively. The velocity U of the collapsing 
fluid is defined as 


U = rDtB < 0 {in the case of collapse). (7) 

Also, the mass function m{t,r) of Cahill and McVittie is obtained from the Riemann tensor component P 23 ^^ 
and is for metric m 


i{t,r) = 


{rB)^ „ „ r3 PP2 r3 B’^ 


p 23 _ 

-A 23 — 


2 A2 2 B 


- - - r^B' 


E is defined as 


E = 


{rB)' 

B 


1 + [/2 _ ‘^rn{t,r) 
rB 


1/2 


Next, the corresponding transport equation for the heat flux reads 


+ < 7 “ = -Kh^^{T ,0 + Tap) - ^ 


( 8 ) 

(9) 


( 10 ) 


where is the projector onto the three space orthogonal to U^, k denotes the thermal conductivity, and T and r 
denote temperature and relaxation time respectively. Observe that due to the symmetry of the problem, equation 
(Cni only has one independent component, which may be written as: 


BDtq 


- — - 2^(1 + 

tE tB t R 


kT 

tE 


[m + LttP] R ^ 


KT^qB 
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Then coupling m to mil one obtains (some misprints in eq.(39) in has been corrected here) 


, \r^rT 7^ /I ^ 7^ EkT' EqB 5 qBEU nET'^qB f t \, 

(/i + T^)(l — a)DtU — Egravi^ — a) + E^yd H—^7^—I -1- ' ' 


2R 


2 At \kT^/’ 


where Egrav and Fhyd are defined by 


( 12 ) 


Egrav = -{n + P)[m + AttPR^] 


(13) 


and 


Fhyd = —E^DrP, 


(14) 


with a given by 


kT 


a = 


r{f^ + P)' 


(15) 


Thus as a tends to 1, the effective inertial mass density of the fluid element tends to zero. Furthermore observe 
that Egrav is also multiplied by the factor (1 — a). Indicating that the effective gravitational attraction on any fluid 
element decreases by the same factor as the effective inertial mass (density). This of course is to be expected, from 
the equivalence principle. It is also worth mentioning that Fhyd is in principle independent (at least explicitly) on 
this factor. 

With these last comments in mind, let us now imagine the following situation. As far as the right hand side of 
is negative, the system keeps collapsing. However, let us assume that the collapsing sphere evolves in such a way that, 
for some region of the sphere, the value of a increases and approaches the critical value of 1. Then, as this process 
goes on, the ensuing decreasing of the gravitational force term would eventually lead to a change of the sign of the 
right hand side of Since that would happen for small values of the effective inertial mass density, that would 

imply a strong bouncing of that part of the sphere, even for a small absolute value of the right hand side of 113- 

In the next section a model will be presented where the effect above appears explicitly. For simplicity we shall 
consider a particular case of the transport equation, corresponding to the so called truncated version, in which case 
the last term on the right of 113 is absent [ij, . In this case, m becomes 


{fj, + P){i — a)DtU — Fgravi^ — a) + Fhyd + 


EkT' 

tB 


EqB 

T 


AqBEU 
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(16) 


III. THE MODEL 

In this section we shall present a numerical model where the decreasing of the effective mass mentioned in the 
previous section will produce a bouncing during the evolution of a dissipative sphere. For simplicity we shall assume 
our fluid to be shear-free and conformally flat, and also that a relevant increase of a takes place only at the boundary 
surface of the sphere. Thus we shall need only to integrate at the boundary surface, implying that we shall deal with 
ordinary differential equations for variables defined on that surface. 


A. The general form of the metric and the field equations 


If the fluid sphere is shear-free and conformally flat, the metric functions take the form 

A= [Ci(t)r2 + 1]H 

and 

P 1 

C'2(t)r2 +C3(t)’ 

where Ci, C 2 and C 3 are arbitrary functions of t. Although the shear free and the conformally flat conditions 
are introduced here in order to simplify calculations, it is worth noticing that these conditions generalize physical 
assumptions widely used in astrophysics. Indeed, the shear free condition in the Newtonian regime describes the 


(17) 

(18) 
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homologous evolution and has been extensively considered in general relativity |^ . On the other hand it is well 
known that the conformally flat condition implies in the perfect fluid case the homogeneity of the energy density 
distribution. 

For the numerical integration we shall need to write all variables in dimensionless form, accordingly we shall redefine 
the metric functions Ci and C 2 by (Ca is already dimensionless): 


where r = = constant defines the boundary surface of the fluid sphere. 

In terms of these dimensionless functions, A and B become 

A=[Ciit)ir/r^)^ + 1]B 

and 

1 

^ C2{t){r/rsy + Csit)' 

Then the following expressions for the physical variables are obtained from Einstein equations 


3 ( C2{r/r^f + C3 
Stt y Ci{r/rj:Y + 1 


+ :^C2C3, 


(19) 


( 20 ) 


( 21 ) 


( 22 ) 


= ^777TT^^^[2(C2(r/rs)2 + C:i){C2{r/r^f + C3) 

8 tt{Ci (r/rs) + 1)"^ 

-3(C2(r/rs)2 + C,f - 2 {C2{r/r^f + C,) 

Ci[r/r^y + 1 

+ ft)| + 2 „(c.(rA„P + l) 

+C3(CiC3 - 2C2)], (23) 


qrl = ^{r/rY.){C3Ci 


( C 2 {r/r^f + C^ V 

V C'i(r/rE)2 + 1 / 


where from now on dot denotes derivative with respect to t/r^. 


(24) 


B. The surface equations 

Next, from © we obtain for the dimensionless proper radius of the sphere 


and from m evaluated at the boundary surface 


Solving these two equations we obtain 


Rn = 


= 


C 2 = 


C.= 


C2 + C3 ’ 

(25) 

C 3 -C 2 

C 3 + C 2 ' 

(26) 

1-Es 

2 i?E 

(27) 

1 + Es 

(28) 

2 i?E 


and 
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On the other hand we have 


As = (Cl + l)i?E, 


(29) 


and from 0 


Us 


C2 + C3 

(Cl + 1 )(C 2 + C 3 ) ’ 


(30) 


where, again, dots denote derivatives with respect to the dimensionless time t/rs- 
may write 


Using 113 and (P|l in pn| we 


Rn = AsUs, 


(31) 


which is our first equation at the surface. 

Next, we use the total loss of mass equation which can be easily derived from (0 and the junction condition 


Ps = {qB)s, 

to obtain (see i) for details) 


(32) 


Ms = —Q{t)AsRs{Us + Es), (33) 

where Ms is the dimensionless mass, Q{t) = AnqsR^^ and qs denotes the dimensionless heat flow qr'^, evaluated at 
the boundary surface. This is our second surface equation. 

Finally, in order to obtain the third surface equation we proceed as follows. From the equations 113 and OH) in 
it can be shown that 


Ce = [Ae (3U| - 1 - - 2RsQ{t)) + 2Es {As - 2Rs)] . (34) 

ZKs 

This is the third equation to be integrated at the surface of the distribution. 

Thus we have a system of three equations (EJ, (j23 and m for the five unknown functions of time Rs,As, Us, Es 
and Q. In order to integrate such a system, we shall prescribe the “luminosity” (Q), and obtain a constraint equation 
from 113, on what we shall elaborate as follows . 

From the dynamic equation 113, using the boundary condition 123 we obtain the pressure gradient at the surface 
Ps = ~ [47ri?E(/rE/3 + Ps) + Us/As] —Ts^, (35) 

where primes and dots denote derivatives with respect to the dimensionless variables r/rs and t/rs respectively, 


Me 


3Me 
47ri?| ’ 


(36) 


Ce 


Q{t) 

AttRs ’ 


(37) 


(observe that ^s and Ps denote the dimensionless expressions for the energy density and pressure evaluated at the 
boundary surface, i.e. these variables multiplied by r|) and 


Ts = as 



Es 

Ry. 


EsQ{t) AQEsUs 

AttRst AttR^ ’ 


(38) 


where we have used kT = ar(/r + P) (conveniently adimensionalised) and have assumed for simplicity that a'^ = 0. 
Finally, using the field equations we obtain the following expression for As 


T Rs/Aasi'iMs / Rs + QRy) — QRy{^ + cke)] 
2ras(l + Es){QRs + iMs/Rs) - QRs{Ry - 7Ust) ' 


( 39 ) 
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C. Strategy of integration 

The integration scheme is now an easy shot: Giving initial conditions for i?s, C/s and My;, and prescribing as and 
Q{t), we can integrate numerically, equations lUTTll . (TTHll . (T^ . with the constraint equation (TTHll . 

The form of as is suggested by the very idea underlying the motivation of this work, namely: the fact that as a 
increases, the ensuing reduction of the gravitational term in the dynamical equation may lead to a bouncing of the 
sphere. Accordingly, we shall take for as a smooth function of time, rising from zero to some value below the critical 
one (a = 1). 


D. Model 

We have ran a large number of models exhibiting bouncing, under physically reasonable conditions, corresponding 
to a wide range of initial data and values of the parameters, and very different choices of Q(t) and a(t). For all these 
choices, the qualitative behaviour associated to the increase of a is essentially the same. From them we have selected 
the following model. 

The initial conditions are 


i?E(0) = 20, 


Ms(0) = 1, 


UsiO) = - 0 . 1 , 

with T = 0.1. These values correspond to a sphere with an initial radius of the order of 400 Km, an initial mass of 
the order of 10 solar masses and a relaxation time of the order of 10“^ seconds. 

The sphere is assumed to be radiating according to 

Q(t) = 

where Qo = 0.001, tm = 0.5 and ct = 0.005, producing a total mass ejection of the order of 0.1%. 

For a we choose 


a(t) = a„/(e-(‘-‘-)/" + 1), 


with am ranging from 0 to 1. 


IV. DISCUSSION 

The influence of pre-relaxation effects on gravitational collapse has been brought out in many works in last decade 
[^ . however the specific effect of bouncing, associated with the decreasing of the effective inertial mass density, 
produced by the increasing of a, had not been illustrated until now. It is worth stressing that a-terms in Eq. m 
come from the inertial factor To/j in Eq. IIOII . 

In this work we provide a numerical model of such bouncing, by assuming an increasing of a at the boundary 
surface. We have concentrated the increase of a on the boundary surface to illustrate the effect, the remaining of the 
sphere is assumed to be dissipating at much lower values of a. Of course, the increasing of a may in principle occur 
at any region of the sphere and even in more that one, simultaneously. The results of our integration is deployed in 
the hguresm^ which exhibit the evolution of different variables with respect to the dimensionless time t/r-^. 

Figure n shows the evolution of i?E for different values of am from 0 to 1, the bouncing is clearly exhibited as well 
as its dependence on a. Figure [21 emphasizes further the link between the increasing of a and the bouncing. 

Figures EHOl shows the behaviour of (dimensionless) energy density, pressure, heat flow and temperature, evaluated 
at the boundary surface. Their values are always regular and satisfy the physical conditions p > P > 0. 

The dimensionless quantity kTs plotted in Figure El is, in conventional units, 

kTs = 2 


( 40 ) 
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with G and c denoting the gravitational constant and the speed of light, and where [k] and [T] denote the numerical 
values of conductivity and temperature in g cm~^ K~^ and K respectively. Therefore the maximum values of 
reached just after the bouncing, correspond to 

(41) 

which may be obtained with [Ts] « 10^^ and [k] « 10^^. These values are well within the acceptable range for those 
variables in a pre-supernovae event |2S|| . 

Thus we have seen that a relatively simple model, whose physical variables exhibit good behaviour and have 
acceptable numerical values, may serve to illustrate the bouncing of a dissipating self-gravitating sphere, produced 
by the decreasing of its effective inertial mass density associated to an increasing of a. 

Nevertheless, in spite of the appeal of the presented model, we are well aware that invoking such an effect to describe 
a specific observed phenomena, would require a much more detailed astrophysical setting. This, however, is out of 
the scope of this paper. 
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FIG. 6: Temperature at the surface (multiplied by kIO^) evolution for different values of am'- 0.0 (solid line); 0.2 (large dashed 
line); 0.4 (short dashed line); 0.6 (dotted line); 0.8 (dot-large dashed line) and 1 (dot-short dashed line). 



